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We  present  a  study  of  internal  pressure  and  surface  tension  of  bare  and  hydrogen  coated  silicon 
nanoparticles  of  2-10  nm  diameter  as  a  function  of  temperature,  using  molecular  dynamics 
simulations  employing  a  reparametrized  Kohen-Tully-Stillinger  interatomic  potential.  The  internal 
pressure  was  found  to  increase  with  decreasing  particle  size  but  the  density  was  found  to  be 
independent  of  the  particle  size.  We  showed  that  for  covalent  bond  structures,  changes  in  surface 
curvature  and  the  associated  surface  forces  were  not  sufficient  to  significantly  change  bond  lengths 
and  angles.  Thus,  the  surface  tension  was  also  found  to  be  independent  of  the  particle  size.  Surface 
tension  was  found  to  decrease  with  increasing  particle  temperature  while  the  internal  pressure  did 
not  vary  with  temperature.  The  presence  of  hydrogen  on  the  surface  of  a  particle  significantly 
reduces  surface  tension  (e.g.,  drops  from  0.83  J/m2  to  0.42  J/m2  at  1500  K).  The  computed  pressure 
of  bare  and  coated  particles  was  found  to  follow  the  classical  Laplace-Young  equation.  ©  2004 
American  Institute  of  Physics.  [DOI:  10.1063/1.1797073] 


I.  INTRODUCTION 

Nanoclusters  of  silicon  are  of  considerable  interest  due 
to  their  potential  application  in  optoelectronics.  At  this  length 
scale,  quantum  confinement  effects  play  an  important  role  in 
the  electronic  properties  which  has  prompted  studies  in  un¬ 
derstanding  the  change  in  electronic  and  transport  properties 
of  materials  and  the  structures  as  a  function  of  size. 1-4  Of 
course,  the  desirability  of  these  materials  places  a  burden  on 
developing  methods  for  preparing  them  in  the  desired  size, 
shape,  and  purity.  Nanoparticles  and  nanocrystals  of  Si  can 
be  prepared  by  plasma  enhanced  chemical  vapor 
deposition. 5-7  Laser  ablation  of  Si  wafer8  and  thermal  vapor¬ 
ization  of  melted  Si  (Ref.  9)  may  also  be  used  to  prepare  Si 
nanoclusters.  In  many  cases  these  methods  result  in  a  silicon 
oxide  layer  on  the  nanoclusters  which  may  be  removed  by 
etching  with  Hydrofluoric  (HF). 

For  production  of  nanoparticles  by  gas-phase  condensa¬ 
tion  processes,  the  high  concentrations  of  fine  particles 
causes  rapid  coagulation.  The  ultimate  size  of  the  spherical 
primary  particles  and  the  growth  of  agglomerates  are  deter¬ 
mined  by  the  competition  between  the  time  for  particle- 
particle  collisions  and  coalescence.10  Particles  will  coalesce 
before  another  collision  event  occurs  when  the  characteristic 
collision  time  between  particles  is  smaller  than  the  charac¬ 
teristic  coalescence  time.  On  the  other  hand,  chain  aggre¬ 
gates  will  be  formed  when  the  collision  time  is  less  than  the 
characteristic  coalescence  time. 

^coalescence^  ^"collision  >  Spherical  particle, 

^coalescence^  ^"collision  ’  Agglomerate. 
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In  principle  then,  by  either  controlling  the  characteristic 
coalescence  time  or  the  collision  time  it  is  possible  to  control 
the  particle  morphology  and  particle  size.  The  two  most  ob¬ 
vious  ways  to  control  the  primary  particle  size  is  to  either 
change  the  characteristic  collision  time  by  dilution  or  change 
the  coalescence  time  by  changing  particle  temperature. 

Another  strategy,  however,  would  be  to  make  particle- 
particle  collisions  nonreactive.  In  a  recent  molecular  dynam¬ 
ics  study  Hawa  and  Zachariah11  showed  that  the  presence  of 
a  hydrogen  passivation  surface  on  silicon  nanoparticles  (size 
200-6400  Si  atoms  at  300-1800  K)  resulted  in  particle 
bounce  during  collision  even  when  the  particles  were  liquid 
droplets.  That  work  also  showed  that  there  was  a  critical 
approach  energy  for  reaction  which  increases  with  increasing 
particle  size  but  decreases  with  increasing  particles  tempera¬ 
ture. 

The  rate  at  which  two  colliding  particles  coalesce  has 
been  modeled  by  both  phenomenological12-17  and  molecular 
dynamics  methods18  which  has  resulted  in  simple  to  apply 
kinetic  laws.  For  example,  the  characteristic  coalescence 
time  calculated  from  a  solid  state  diffusion  model  is  written 
as12’13'15 

3  kTpN 

647 tctD'  ^ 

where  Tp  is  the  particle  temperature,  N  is  the  number  con¬ 
centration,  D  is  the  solid-state  diffusion  coefficient  reported 
as  an  Arrhenius  function  of  the  temperature,16  and  cr  is  the 
surface  tension  of  the  particle.  For  viscous  flow,  the  coales¬ 
cence  time  is  given  by17 
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where  dp  is  the  diameter  of  the  particle  and  77  is  the  tempera¬ 
ture  dependent  viscosity. 

In  looking  at  the  above  two  kinetic  laws  it  is  apparent 
that  surface  properties  play  an  important  role.  The  driving 
force  for  coalescence  is  in  fact  the  minimization  of  the  sur¬ 
face  free  energy  and  it  is  not  surprising  then  that  surface 
tension  appears  as  a  parameter  in  both  equations,  and  the 
diffusion  constant  is  obviously  important  in  the  solid-state 
coalescence  event. 

Given  that  the  driving  force  for  coalescence  is  the  mini¬ 
mization  of  the  surface  free  energy,  it  is  important  to  under¬ 
stand  the  effect  of  hydrogen  surface  passivation  on  the  sur¬ 
face  tension  and  the  diffusion  coefficient  of  the  particles  as 
well  as  the  role  of  particle  size.  Although  hydrogenated 
amorphous  silicon  has  found  a  variety  of  technological  ap¬ 
plications,  and  has  been  studied  extensively,19  these  previous 
studies  focused  on  the  electronic  and  optical  properties, 20-26 
thin  film  growth,27'28  and  the  deposition  of  clusters  on  hydro¬ 
genated  surface.29-38  To  our  best  knowledge,  surface  tension 
and  diffusion  coefficient  of  hydrogenated  Si  nanoparticles 
have  not  been  reported  on. 

In  this  paper  we  investigate  the  role  of  hydrogen  passi¬ 
vation  on  physical  properties,  such  as  surface  tension,  inter¬ 
nal  pressure,  and  surface  diffusion.  We  use  classical  molecu¬ 
lar  dynamics  simulation  using  the  reparametrized  Kohen- 
Tully-Stillinger  (KTS)  potential  for  the  silicon-hydrogen 

1 1  IQ 

system.  ’ 

II.  COMPUTATIONAL  MODEL 
AND  SIMULATION  PROCEDURE 

Recent  geometric  optimizations  and  molecular  dynamics 
(MD)  simulations  based  on  an  empirical  tight-binding  ap¬ 
proach  for  fully  and  partially  hydrogenated  Si  nanocrystals 
showed  that  the  structural  properties  were  not  very  size  de¬ 
pendent,  unlike  electronic  properties.40-41  However,  hydroge¬ 
nation  has  been  shown  to  stabilize  the  structure  of  Si„  clus¬ 
ters  (n  =  2  to  10)  through  termination  of  dangling  bonds.42 
Prior  studies  using  molecular  dynamics  have  been  used  to 
study  surface  energies  of  nanoclusters.  For  example,  studies 
of  solid  and  liquid  silica  particles  show  that  the  surface  ten¬ 
sion  was  independent  of  particle  size  over  a  wide  range  of 
temperature.43  On  the  other  hand,  MD  studies  of  Lennard- 
Jones  (LJ)  clusters44  and  a  variable  charge  transfer  MD  for 
aluminum45  showed  that  surface  tension  decreased  with  de¬ 
creasing  particle  size. 

This  study  involves  atomistic  simulations  using  classical 
molecular  dynamics.  For  this  study  we  use  the  reparam¬ 
etrized  KTS  interatomic  potential  for  the  silicon-hydrogen 
system  developed  by  Hawa  and  Zachariah. 11  This  inter¬ 
atomic  potential  for  silicon  was  originally  developed  by 
Stillinger  and  Weber46  and  extended  by  Kohen  et  al.39  to 
include  Si-H  and  H-H  interactions.  Similar  sets  of  potential 
energy  functions  have  also  been  developed  by  Murty  and 
Atwater,  Ohira  and  co-workers  ’  and  Ramalingam 
et  al.31  where  a  Tersoff-type  potential48-51  was  extended  to 
describe  interatomic  interactions  in  the  Si:H  system.  This 
extended  version  of  the  Tersoff  potential  has  been  tested  suc¬ 
cessfully  for  its  accuracy  in  describing  the  Si:H  system  in 
several  earlier  studies;  however,  the  simulation  of  liquid  sili¬ 


con  has  not  been  well  described  by  the  potential.48  By  con¬ 
trast,  the  reparametrized  KTS  potential,  which  is  the  ex¬ 
tended  spherical  wave  Stillinger  and  Weber  (SW)  potential, 
was  designed  to  describe  interactions  in  both  solid  and  liquid 
forms  of  silicon.  Since  most  synthesis  processes  leading  to 
cluster  formation  occur  at  high  temperature,  cluster  growth 
by  coalescence  is  dominated  by  liquidlike  characteristics, 
and  the  accuracy  of  the  SW  potential  increases  with  increas¬ 
ing  particle  size  or  temperature,  we  use  this  potential  form 
for  our  investigations.  The  reparametrized  KTS  potential  en¬ 
ergy  is  a  sum  of  two  and  a  three -body  interactions,  and  the 
details  of  the  model  and  its  parameters  are  given  in  the  lit¬ 
erature  (Hawa  and  Zachariah).11 

All  simulations  were  run  on  an  Origin  computer  running 
up  to  16  processors.  The  trajectories  of  all  the  atoms  are 
determined  by  integrating  the  equation  of  motion  according 
to  the  velocity  form  of  the  Verlet  algorithm52  with  rescaling 
of  atomic  velocities  at  each  time  step  to  achieve  temperature 
control.  Time  steps  of  0.5  and  0.05  fs  were  typically  used  for 
both  pure  silicon,  and  for  hydrogen  coated  particles  to  ensure 
energy  conservation,  and  the  Verlet  neighbor  list  with  paral¬ 
lel  architecture  was  employed  in  all  the  simulations.  The 
neighbor  list  was  renewed  every  ten  steps.  The  simulations 
take  place  in  a  spherical  cavity  of  10  nm  radius  with  an 
elastic  boundary  condition. 

The  first  step  in  the  equilibration  process  was  to  prepare 
pure  silicon  particles  of  various  sizes  (2,  3,  4,  6,  and  10  nm 
diameters)  (200,  800,  1600,  6400,  25  600  Si  atoms)  at  300  K. 
After  the  angular  momentum  was  removed,  the  particle  tem¬ 
perature  was  raised  to  2100  K  using  constant  temperature 
MD  for  50  ps.  The  particle  temperature  was  then  reduced 
slowly  to  300  K  and  equilibrated  for  50  ps.  The  next  step 
was  to  coat  the  particles  with  hydrogen  atoms.  Since  the 
particles  were  already  equilibrated,  almost  all  surface  atoms 
had  a  coordination  number  of  three.  A  hydrogen  atom  was 
placed  on  each  surface  silicon  atoms  at  a  distance  of  1 .5  A  to 
prepare  fully  hydrogen  coated  particles,  while  the  particle 
temperature  was  maintained  at  300  K  for  10  ps.  Any  hydro¬ 
gen  atoms  that  were  released  from  the  surface  were  removed 
from  the  simulation,  and  the  dynamics  repeated  for  10  ps.  As 
a  result,  the  numbers  of  hydrogens  placed  on  the  silicon  par¬ 
ticles  composed  of  200,  800,  1600,  and  6400  Si  atoms  was 
74,  211,  367,  and  785,  respectively.  For  the  preparation  of 
partially  hydrogen  coated  particles,  appropriate  number  of 
hydrogen  atoms  were  randomly  removed  from  the  particles, 
and  the  dynamics  repeated  for  10  ps  for  equilibration.  After 
generating  the  hydrogen  monolayer  on  the  silicon  particles, 
the  temperature  of  the  particles  was  slowly  raised  to  the  de¬ 
sired  temperatures  of  1000,  1500,  or  1800  K  and  maintained 
at  constant  temperature  for  50  ps.  For  the  last  step  in  the 
preparation  process,  the  simulations  were  switched  to  a  con¬ 
stant  energy  calculation  for  20  ps.  If  the  average  temperature 
of  the  particle  deviated  by  more  than  10  K  over  this  period, 
the  equilibration  process  was  repeated  until  the  particle  tem¬ 
perature  deviated  by  less  than  10  K.  An  external  and  cross- 
sectional  view  of  an  equilibrated  nanoparticle,  consisting  of 
6400  Si  atoms,  coated  with  785  H  atoms,  corresponding  to 
roughly  a  6  nm  particle  are  shown  in  Fig.  1.  The  properties 
such  as  radial  pressure,  density,  surface  tension,  and  diffusiv- 
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FIG.  1.  A  6  nm  hydrogen-coated  silicon  nanoparticles  (6400Si+785  H 
atoms):  (a)  external  view,  (b)  cross-sectional  view. 

ity  were  computed,  and  the  pressure  and  density  were  aver¬ 
aged  by  collecting  200  snapshots  over  40  ps. 

III.  RESULTS  AND  DISCUSSION 
A.  Density 

The  radial  density  distributions  for  fully  and  50%  hydro¬ 
gen  coated  and  bare  4  nm  silicon  particles  at  1500  K,  is 
shown  in  Fig.  2.  The  hydrogen  density  profiles  of  both  100% 
and  50%  coated  particles  show  that  the  maximum  values  of 
the  hydrogen  density  is  achieved  at  the  surface  of  the  silicon 
particle  and  that  no  hydrogen  exists  inside  of  either  particle. 
The  peak  value  in  hydrogen  densities  occurs  at  the  same 
radial  location  independent  on  hydrogen  coverage.  Naturally 
the  magnitude  in  density  of  the  50%  coating  is  lower  than 
that  of  a  fully  coated  particle  and  the  ratio  of  the  areas  under 
the  distribution  curves  of  the  fully  coated  particle  to  that  of  a 
50%  coated  particle  is  two.  The  fact  that  the  higher  coverage 
condition  seemed  to  have  a  density  profile  that  was  wider, 
and  biased  toward  the  interior  of  the  particle,  might  indicate 
that  the  silicon  atoms  that  have  hydrogen  attached  are  more 
likely  to  sit  out  of  the  surface,  while  the  unpassivated  atoms 
recede  toward  the  interior.  On  the  other  hand,  the  silicon 
density  profiles  for  all  bare,  fully,  and  partially  coated  par¬ 
ticles  are  indistinguishable  except  near  the  surface.  Prior 
work  using  an  empirical  tight-binding  approach  for  studying 
the  structure  of  fully  and  partially  hydrogenated  Si  nanocrys¬ 
tals  (Si<706  atoms)  found  that  there  was  only  a  very  small 
lattice  strain,  of  the  order  1 0  4  to  10  3 , 4 1  in  agreement  with 
x-ray  diffraction  measurements,40  and  confirms  that  the  den¬ 
sity  profiles  should  be  independent  of  surface  coverage.  Fig- 


FIG.  2.  Radial  density  profiles  of  fully  and  50%  hydrogen  coated  and  bare 
4  nm  silicon  particle  at  1500  K. 
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FIG.  3.  Radial  density  profile  for  a  3  and  6  nm  hydrogen  coated  silicon 
particle  at  1800  K. 

ure  3  shows  the  density  profiles  of  3  and  6  nm  bare  silicon 
particles  at  1800  K.  It  can  be  seen  that  the  average  densities 
do  not  depend  on  the  size  of  the  particles.  The  average  den¬ 
sity  of  the  particle  at  1800  K  is  about  2.57  g/cc,  and  the  value 
is  consistent  with  the  density  recently  measured  for  bulk  pure 
molten  silicon  (2.55  g/cc).53  Our  studies  show  neither  par¬ 
ticle  size  or  hydrogen  surface  coverage  have  an  effect  on 
silicon  particle  density. 

B.  Diffusion 

Self-diffusion  coefficient  D  is  calculated  from  the  slope 
of  the  atomic  mean-square  displacement  (MSD): 

1  N 

MSD=-2  [o(f)-o(0)]2,  (3) 

i=  J 

D=^<MSCh  (4) 

where  N  is  the  number  of  atoms  in  the  particle  and  r  is  the 
position  of  each  atom.  Since  the  curvature  in  the  MSD  to¬ 
ward  an  asymptote  results  from  the  finite  nature  of  the  clus¬ 
ter,  the  self-diffusion  coefficient  can  be  evaluated  from  the 
initial  slope  of  the  MSD  in  Eq.  (3).  Figure  4  shows  the  self¬ 
diffusion  coefficients  of  both  bare  and  fully  hydrogenated  6 
nm  particles  at  1000,  1500,  and  1800  K.  We  find  that  the 
presence  of  hydrogen  on  the  surface  of  a  silicon  particle 
slightly  increases  the  diffusivity  of  silicon  atoms  of  a  coated 


FIG.  4.  Self-diffusion  coefficients  of  both  bare  and  coated  6  nm  particles  at 
1000,  1500,  and  1800  K. 
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particle  by  about  5-15%  over  the  temperature  range 
sampled.  Since  our  reparametrized  KTS  potential  is  an  ex¬ 
tension  of  SW  potential,  our  results  show  agreement  with  the 
earlier  work  by  Zachariah  el  al ,54  who  used  the  SW  potential 
for  Si  clusters  up  to  480  atoms.  By  comparing  the  simulation 
results  with  the  experimental  data  of  surface  self¬ 
diffusion,55  57  they  found  substantial  discrepancy  in  the  bulk 
experimental  data.  The  SW  potential  is  known  to  underpre¬ 
dict  the  binding  energy  in  bulk  silicon  that  would  tend  to 
increase  the  computed  surface  diffusion  coefficient.  How¬ 
ever,  their  results  showed  agreement  with  those  of  tight- 
binding  molecular  dynamics.58 

The  self-diffusion  coefficient  of  hydrogen  atoms  is  also 
shown  in  Fig.  4.  We  find  that  hydrogen  diffusivity  is  essen¬ 
tially  the  same  as  that  of  silicon  diffusivity.  Variable  tem¬ 
perature  scanning  tunnel  microscope  (STM)  measurements 
on  low  hydrogen  coverage  silicon  films  at  temperatures  rang¬ 
ing  from  300  to  700  K  showed  that  individual  hydrogen 
atoms  became  mobile  at  around  570  K,  and  that  rate  of  hop¬ 
ping  along  the  dimer  rows  was  consistent  with  an  activation 
energy  of  1.68  ±  0. 15  eV.59  In  a  related  study,  an  upper 
bound  of  D(7’=740K)  <l(U9cm2/s  was  obtained.60 
Reider  et  al.  employed  an  altogether  different  experimental 
approach  and  found  effectively  the  same  barrier  height  ( Ea 
=  1.5±0.2eV)  and  a  preexponential  factor  D0  of 
1 0  3  cm2/s  using  the  method  of  second-harmonic  diffraction 
from  adsorbate  gratings.61  Using  these  experimental  values 
we  estimate  a  diffusion  coefficient  (hopping)  for  H  at  1000  K 
of  ~2X  10  11  cm2/s,  which  is  much  lower  than  our  simu¬ 
lation  results  of  the  self-diffusion  coefficient  of  hydrogen 
atoms  in  Fig.  4.  This  implies  that  the  surface  hydrogen  mo¬ 
bility  directly  reflects  the  mobility  of  its  silicon  anchor,  and 
H  hopping  rate  is  negligible  at  our  temperature  range,  rela¬ 
tive  to  Si  self-diffusion. 

C.  Pressure 

The  internal  pressures  are  computed  using  the  Irving- 
Kirkwood  pressure  tensor  extended  to  a  spherical  symmetric 
system  as  shown  by  Thompson.44  The  normal  component  to 
the  spherical  surface  of  the  pressure  Ppj(r)  is  given  by 

PN(r)  =  PK(r)  +  Pu(r),  (5) 

where  Px(r)  and  Pu(r)  are  kinetic  and  configurational 
terms,  respectively,  and  are  given  by 

P  K{r)  =  kBTp(r), 

(6) 

Pu(r)  =  S~1'Z  fk, 

k 

where  kB  is  Boltzmann  constant,  T  is  a  particle  temperature, 
S  is  the  area  of  spherical  surface  of  radius  r,  and  fk  is  the 
normal  component  of  the  force  between  an  atomic  pair  act¬ 
ing  across  the  surface  S.  We  introduce  800  subspherical  cal¬ 
culation  surfaces  in  the  spherical  computational  cavity  of  10 
nm  radius,  to  compute  the  configurational  term  of  pressure 
for  every  subsphere  surface.  The  accuracy  of  this  method  of 
pressure  computation  depends  on  the  number  of  atoms 
within  the  subsphere,  and  therefore  quickly  diminishes  with 
decreasing  radius. 


Radius  (A) 

FIG.  5.  Radial  pressure  profile  for  2-10  nm  silicon  particles  at  1500  K. 


The  radial  internal  pressure  distributions  for  2  to  10  nm 
particles  at  1500  K  are  shown  in  Fig.  5.  We  see  that  the 
internal  pressure  distributions  are  very  sensitive  to  the  size  of 
the  particle  at  a  fixed  temperature.  The  internal  pressure  in¬ 
creases  with  decreasing  particle  size,  consistent  with  the  in¬ 
ternal  pressure  profile  reported  for  silica.43  The  effect  of  tem¬ 
perature  are  reported  in  Fig.  6,  and  show  that  internal 
pressure  is  independent  of  the  particle  temperature,  even 
though  the  kinetic  contribution  to  the  normal  pressure  PN 
depends  linearly  on  the  temperature.  This  result  is  also  con¬ 
sistent  with  the  study  for  silica  reported  in  the  literature.43 
The  details  of  these  results  and  its  relationship  to  the  well 
known  Laplace-Young  equation  will  be  differed,  so  as  to 
present  results  on  surface  tension  as  these  effects  are  related. 
The  radial  pressure  distribution  at  1000  K  shows  oscillations 
near  the  core  of  the  particle.  This  is  a  manifestation  of  the 
fact  that  the  number  of  atoms  within  the  subsphere  calcula¬ 
tion  surface  decreases  as  the  inverse  cube  of  the  radial  dis¬ 
tance  inward,  thus  providing  poor  statistics,  and  the  fact  that 
since  the  particles  are  solidlike  the  atoms  have  poor  mobility. 
If  we  could  integrate  for  long  enough,  one  should  expect  to 
see  the  oscillations  disappear  as  they  have  for  liquidlike  par¬ 
ticles. 

Figure  7  shows  the  pressure  profiles  of  bare,  50%,  and 
fully  coated  6  nm  particles  at  1500  K.  The  presence  of  hy¬ 
drogen  on  the  surface  reduces  the  internal  pressure  of  the 
coated  particle  by  about  40%.  In  the  definition  of  the 
Irving-Kirkwood  pressure  tensor,  the  attractive  forces  make  a 


FIG.  6.  Radial  pressure  profile  for  6  nm  silicon  particles  at  1000  and 
1500  K. 


Downloaded  03  Nov  2004  to  129.2.63.242.  Redistribution  subject  to  AIP  license  or  copyright,  see  http://jcp.aip.org/jcp/copyright.jsp 


9047 


J.  Chem.  Phys.,  Vol.  121,  No.  18,  8  November  2004 


Internal  pressure  and  surface  tension  of  nanoparticles 


FIG.  7.  Radial  pressure  profile  for  6  nm  bare  and  50%  and  fully  coated 
silicon  particles  at  1500  K. 


FIG.  9.  Surface  tension  of  bare  silicon  particles  as  a  function  of  particle  size 
at  1000,  1500,  and  1800  K. 


negative  contribution  to  the  pressure,  while  repulsive  forces 
make  a  positive  contribution.  Thus  a  decrease  in  pressure 
indicates  that  the  bond  lengths  of  all  atoms  in  the  particle  are 
on  average  elongated  by  the  presence  of  hydrogen  on  the 
surface  of  the  particles. 

Figure  8  shows  that  the  average  bond  length  for  surface 
atoms  for  both  bare  and  coated  silicon  as  a  function  of  par¬ 
ticle  size  at  1500  K.  We  can  clearly  see  that  the  average  bond 
lengths  of  hydrogen  coated  particles  are  longer  than  those  of 
bare  particles  for  all  sizes.  Moreover,  both  average  surface 
bond  length  initially  increases  with  increasing  particle  size, 
but  that  the  size  dependence  disappears  by  about  4  nm 
(—2000  atoms).  It  is  known  that  the  average  bond  length 
tends  to  become  shorter  with  decrease  of  the  coordination 
number.62  The  similar  tendency  in  the  average  coordination 
number  has  been  also  pointed  out  in  our  previous  studies 
(Zachariah  et  al.54). 

The  average  bond  lengths  of  inner  1800  silicon  atoms  of 
bare  and  fully  coated  particles  at  1500  K  are  2.4928  and 
2.4944  A,  and  the  ratio  of  the  difference  to  the  bond  length 
of  bare  particles  is  6.4X  10  4,  in  agreement  with  the  strain 
of  1 0  4  to  1 0  3  obtained  from  tight-binding  calculation  of 
silicon  crystals41  and  x-ray  diffraction  measurements.40 

D.  Surface  tension 

The  surface  tension  cr  is  computed  using  the  mechanical 
approach  reported  by  Thompson  et  al.44 


FIG.  8.  Average  bond  length  of  bare  and  hydrogen  coated  surface  silicon 
atoms  as  a  function  of  particle  size  at  1500  K. 
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where  pc  is  gas  pressure  outside  the  particle  and  pL  is  pres¬ 
sure  inside  the  particle.  Since  almost  all  contribution  to  the 
surface  tension  comes  from  the  surface  pressure,  we  esti¬ 
mated  the  pL  as  the  average  pressure  near  the  surface  of  the 
particle  and  pc  can  be  ignored. 

The  surface  tension  as  a  function  of  particle  size  is 
shown  in  Fig.  9.  It  is  interesting  to  note  that  the  surface 
tension  does  not  depend  on  particle  size,  and  is  qualitatively 
consistent  with  the  study  for  silica.43  Since  the  surface  ten¬ 
sion  depends  on  both  the  calculated  pressure  and  density 
profiles,  the  accuracy  of  the  surface  tension  diminishes  with 
decreasing  particle  size. 

Thompson  et  al.44  and  Sonwane  et  al.45  in  their  work 
reported  that  surface  tension  for  Lennard-Iones  clusters  and 
aluminum,  respectively,  decreased  with  decreasing  particle 
size.  We  believe  the  distinction  between  LJ  and  aluminum 
versus  silicon  and  silica  can  be  attributed  to  the  nature  of 
bonding.  Silica  and  silicon  are  covalent  structures,  whose 
bonding  is  highly  directional.  Changes  in  surface  curvature 
and  the  associated  surface  forces  are  not  sufficient  to  signifi¬ 
cantly  change  bond  lengths  and  angles.  This  implies  that  the 
nature  of  covalent  bonded  surfaces  is  to  first  order  size  in¬ 
variant.  In  contrast  LJ  and  possibly  metals  are  more  able  to 
find  different  bonding  configurations  as  the  surface  curvature 
changes,  thus  resulting  in  changes  to  surface  energy  as  a 
function  of  particle  size. 

In  order  to  evaluate  this  argument,  we  ignore  the  third- 
body  effect  in  the  reparametrized  KTS  potential  so  that  the 
potential  behaves  like  an  LJ  cluster.  Figures  10  and  11  show 
the  radial  internal  pressure  distribution  for  various  particle 
sizes,  and  the  surface  tension  as  a  function  of  the  particle  size 
at  fixed  temperature.  Both  the  pressure  and  surface  tension 
values  are  normalized  by  the  maximum  values  obtained  in 
the  simulation.  The  internal  pressures  obtained  are  qualita¬ 
tively  similar  to  that  obtained  with  the  full  potential  (see  Fig. 
5).  However,  the  surface  tension  shows  a  size  dependence 
that  is  similar  to  that  observed  by  Thompson  et  al.  for  LJ 
clusters  and  Sonwane  et  al.  for  aluminum.  This  implies  that 
indeed,  consistent  with  our  conjecture,  directionally  bonded 
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FIG.  10.  Radial  pressure  profile  neglecting  the  third  body  term  for  various 
sizes  for  at  1500  K. 

materials  should  not  show  any  size  dependent  surface  ten¬ 
sion  effects. 

From  Fig.  9  we  also  see  that  the  surface  tension  de¬ 
creases  with  increasing  particle  temperature  for  the  pure  sili¬ 
con  particles.  This  trend  is  also  observed  for  hydrogenated 
particles  shown  in  Fig.  12.  Since  we  have  shown  that  surface 
tension  is  size  invariant,  results  presented  in  Fig.  12  are  the 
averages  of  results  for  2  to  10  nm  particles.  The  surface 
tension  is  the  energy  required  to  incrementally  increase  the 
surface  area.  Since  it  is  easier  to  increase  surface  area  at 
higher  temperatures  when  the  atoms  are  more  mobile,  the 
surface  tension  must  decrease  with  increasing  temperature. 
Experimental  results63  of  surface  tension  of  bulk  silicon  as  a 
function  of  temperature  are  also  plotted  in  Fig.  12  and  show 
the  same  temperature  dependence  and  are  of  similar  magni¬ 
tude  (-15%). 

Figures  9  and  12  also  show  the  role  of  hydrogen  atoms 
on  surface  tension.  With  increasing  hydrogen  coverage  the 
surface  tension  decreases.  Conceptually  the  surface  tension  is 
a  bulk  measure  of  the  energy  of  a  surface  atom  relative  to 
interior  atoms.  At  the  atomic  scale  the  attractive  forces  be¬ 
tween  atoms  interior  to  a  particle  are  shared  with  all  neigh¬ 
boring  atoms.  However,  those  on  the  surface  have  fewer 
neighboring  atoms  above  the  surface,  and  exhibit  stronger 
attractive  forces  upon  their  nearest  surface  and  subsurface 
neighbors.  Since  this  enhancement  of  the  interatomic  attrac¬ 


FIG.  11.  Surface  tension  neglecting  the  third  body  term  for  various  sized 
particles  at  1500  K. 
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FIG.  12.  Surface  tension  of  bare  and  partially  and  fully  coated  silicon  par¬ 
ticles  with  experimental  results  and  the  case  of  a  half  Si-H  bonding  energy. 

tive  forces  at  the  surface  is  called  as  surface  tension,  the 
presence  of  hydrogen  atoms  on  the  silicon  surface  produces 
attractive  forces  above  the  silicon  surface  and  reduces  the 
surface  tension  of  the  silicon  particle. 

Since  the  surface  tension  is  a  measure  of  the  energy  of 
surface  atoms,  which  can  be  changed  by  the  nature  of  the 
surface  termination  bond,  one  can  in  principle  tune  surface 
tension  through  attachment  of  an  appropriate  ligand.  For  ex¬ 
ample,  if  we  arbitrarily  change  the  bond  energy  of  Si-H  to 
half  its  nominal  value,  we  get  the  result  shown  by  the  symbol 
(X)  in  Fig.  12,  which  is  an  increase  in  surface  tension  of 
about  30%.  Since  the  surface  tension  is  a  key  parameter  in 
the  kinetics  of  nanoparticle  coalescence  as  described  by  Eqs. 
(1)  and  (2),  the  use  of  surface  manipulation  may  afford  a  new 
opportunity  which  to  this  point  has  been  unexplored  to  affect 
growth  rates  of  gas  phase  generated  nanoparticles. 

Finally  we  turn  our  attention  to  the  well  known 
Laplace-Young  equation,  which  is  widely  used  for  droplets 
to  describe  the  relationship  between  pressure,  surface  ten¬ 
sion,  and  size: 

FL-Fo  =  2o-/r.  (8) 

Figure  13  shows  the  pressure  change  as  a  function  of  radius 
using  Eq.  (8)  at  1500  K  and  compares  it  with  the  internal 
pressure  obtained  from  Fig.  7.  We  can  see  that  the  internal 
pressure  computed  from  Laplace-Young  equation  is  in  ex¬ 
cel- 


FIG.  13.  Comparison  of  pressure  obtained  from  Laplace- Young  equation 
with  Irving-Kirkwood  pressure  tensor. 
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lent  agreement  with  our  MD  results  for  both  the  coated  and 
uncoated  cases,  with  the  exception  of  very  small  particle 
sizes.  The  results  from  Fig.  13  also  imply  the  Laplace- 
Young  equation  is  valid  for  surface  coated  particles,  i.e.,  non- 
homogeneous  particles.  Extrapolation  of  Eq.  (8)  to  very 
small  clusters  imply  increasing  internal  pressure  to  the  point 
of  a  pressure  singularity.  In  fact,  the  concept  of  a  surface 
tension  would  itself  come  into  question  for  very  small  par¬ 
ticle  sizes.  For  the  purposes  of  this  study  we  conclude  that, 
over  the  range  of  2-10  nm,  the  phenomenological  descrip¬ 
tion  provided  by  the  Laplace-Young  equation  is  valid. 

IV.  CONCLUSION 

Classical  molecular  dynamics  using  the  reparametrized 
three-body  KTS  potential  was  used  to  study  the  internal  pres¬ 
sure  and  surface  tension  of  bare  and  hydrogen  coated  silicon 
nanoparticles.  Simulations  were  performed  over  a  wide 
range  of  particle  diameters  between  2-10  nm  and  at  a  tem¬ 
perature  of  1000,  1500,  and  1800  K.  It  was  found  that  inter¬ 
nal  pressure  increases  with  decreasing  particle  size  but  is 
independent  of  particle  temperature.  Surface  tension  de¬ 
creases  with  increasing  particle  temperature  but  is  indepen¬ 
dent  of  the  particle  size.  We  believe  that  this  result  should  be 
generic  to  covalent  bonded  structures,  whose  bond  lengths 
and  angles  are  not  easily  changed,  by  changes  in  surface 
curvature  and  associated  surface  forces.  It  was  observed  that 
the  presence  of  hydrogen  on  the  silicon  particle  has  little 
effect  on  surface  diffusion,  but  significantly  decreases  both 
internal  pressure  as  well  as  surface  tension.  The  internal 
pressure  of  both  bare  and  coated  particles  was  found  to  agree 
with  Laplace-Young  equation. 
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